library(package = "tidyverse")
## ── Attaching packages ───────────────────────────────────────────────────────────────────────────────────── tidyverse 1.2.1 ──
## ✔ ggplot2 3.1.0 ✔ purrr 0.3.2
## ✔ tibble 2.1.1 ✔ dplyr 0.8.0.1
## ✔ tidyr 0.8.3 ✔ stringr 1.4.0
## ✔ readr 1.3.1 ✔ forcats 0.4.0
## ── Conflicts ──────────────────────────────────────────────────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
library(package = "lubridate")
##
## Attaching package: 'lubridate'
## The following object is masked from 'package:base':
##
## date
library(package = "rlist")
library(package = "ncdf4")
library(package = "ncdf4.helpers")
library(package = "PCICt")
library(package = "reshape2")
##
## Attaching package: 'reshape2'
## The following object is masked from 'package:tidyr':
##
## smiths
# directory/URL root
directory = "/projects/BIG_WEATHER/GENS_ERROR_RT/triangle_archives/"
directory = "http://kyrill.ias.sdsmt.edu:8080/thredds/fileServer/BWW_GENS/CI_STAT/"
time_range = "2015-06-01_00_to_2019-03-31_00"
region = "WRFRAP"
RdataFile = str_c(directory,
"gens_03_ensemble__",
"T2M_MSLP_M10_SPCH2M_ISOHGT_U10_V10_FRICV_GUST",
"__error__",
region,
"__",
time_range,
".RData",
sep = "")
load(file = url(RdataFile))
Time = unique(CI_Ensemble_Stats$Time)
Fx_Hour = unique(CI_Ensemble_Stats$Fx_Hour)
Variable = unique(CI_Ensemble_Stats$Variable)
Height = unique(CI_Ensemble_Stats$Height)
CI_Ensemble_Stats$Month = month(CI_Ensemble_Stats$Time)
CI_Ensemble_Stats = CI_Ensemble_Stats %>%
mutate(Quarter = case_when(((Month == 12) | (Month == 01) | (Month == 02))
~ "DJF",
((Month == 03) | (Month == 04) | (Month == 05))
~ "MAM",
((Month == 06) | (Month == 07) | (Month == 08))
~ "JJA",
((Month == 09) | (Month == 10) | (Month == 11))
~ "SON") )
Var = "T2M"
Varname = "2-m Air Temperature"
Hgt = 2
Fx = 24
for (Fx in Fx_Hour[2:length(Fx_Hour)]) {
fx_subset = CI_Ensemble_Stats %>% filter((Variable == Var) &
(Height == Hgt) &
(Fx_Hour == Fx) ) %>%
select(-Quarter)
seasonal_subset = CI_Ensemble_Stats %>% filter((Variable == Var) &
(Height == Hgt) &
(Fx_Hour == Fx) )
myplot = ggplot(data = seasonal_subset) +
aes(x = Ens_StDev,
y = RMSE_Ens000,
color = Quarter) +
facet_wrap(facets = ~ Quarter) +
theme_bw() +
theme(strip.background = element_rect(fill=NA),
) +
labs(title = str_c(Fx,
"-hr Forecast, CI Seasonal Triangles for ",
Varname,
sep = ""),
subtitle = str_c(region_name,
sep = "")) +
xlab("Ensemble Standard Deviation") +
ylab("Root Mean Squared Error") +
geom_point(data = fx_subset,
color = "grey",
alpha = 0.5) +
scale_colour_manual(values=c("DJF" = "cyan",
"MAM" = "green",
"JJA" = "magenta",
"SON" = "orange"),
guide =FALSE) +
geom_point(data = seasonal_subset,
alpha = 0.7)
print(myplot)
}
## Warning: Removed 332 rows containing missing values (geom_point).
## Warning: Removed 83 rows containing missing values (geom_point).
## Warning: Removed 552 rows containing missing values (geom_point).
## Warning: Removed 138 rows containing missing values (geom_point).
## Warning: Removed 748 rows containing missing values (geom_point).
## Warning: Removed 187 rows containing missing values (geom_point).
## Warning: Removed 904 rows containing missing values (geom_point).
## Warning: Removed 226 rows containing missing values (geom_point).
## Warning: Removed 1128 rows containing missing values (geom_point).
## Warning: Removed 282 rows containing missing values (geom_point).
## Warning: Removed 1312 rows containing missing values (geom_point).
## Warning: Removed 328 rows containing missing values (geom_point).
## Warning: Removed 1484 rows containing missing values (geom_point).
## Warning: Removed 371 rows containing missing values (geom_point).
## Warning: Removed 1620 rows containing missing values (geom_point).
## Warning: Removed 405 rows containing missing values (geom_point).
## Warning: Removed 1844 rows containing missing values (geom_point).
## Warning: Removed 461 rows containing missing values (geom_point).
## Warning: Removed 2008 rows containing missing values (geom_point).
## Warning: Removed 502 rows containing missing values (geom_point).
## Warning: Removed 2196 rows containing missing values (geom_point).
## Warning: Removed 549 rows containing missing values (geom_point).
## Warning: Removed 2332 rows containing missing values (geom_point).
## Warning: Removed 583 rows containing missing values (geom_point).
## Warning: Removed 2528 rows containing missing values (geom_point).
## Warning: Removed 632 rows containing missing values (geom_point).
## Warning: Removed 2696 rows containing missing values (geom_point).
## Warning: Removed 674 rows containing missing values (geom_point).
Var = "ISOHGT"
Varname = "500-hPa Isobaric Heights"
Hgt = 500 * 100
Fx = 24
for (Fx in Fx_Hour[2:length(Fx_Hour)]) {
fx_subset = CI_Ensemble_Stats %>% filter((Variable == Var) &
(Height == Hgt) &
(Fx_Hour == Fx) ) %>%
select(-Quarter)
seasonal_subset = CI_Ensemble_Stats %>% filter((Variable == Var) &
(Height == Hgt) &
(Fx_Hour == Fx) )
myplot = ggplot(data = seasonal_subset) +
aes(x = Ens_StDev,
y = RMSE_Ens000,
color = Quarter) +
facet_wrap(facets = ~ Quarter) +
theme_bw() +
theme(strip.background = element_rect(fill=NA),
) +
labs(title = str_c(Fx,
"-hr Forecast, CI Seasonal Triangles for ",
Varname,
sep = ""),
subtitle = str_c(region_name,
sep = "")) +
xlab("Ensemble Standard Deviation") +
ylab("Root Mean Squared Error") +
geom_point(data = fx_subset,
color = "grey",
alpha = 0.5) +
scale_colour_manual(values=c("DJF" = "cyan",
"MAM" = "green",
"JJA" = "magenta",
"SON" = "orange"),
guide =FALSE) +
geom_point(data = seasonal_subset,
alpha = 0.7)
print(myplot)
}
## Warning: Removed 340 rows containing missing values (geom_point).
## Warning: Removed 85 rows containing missing values (geom_point).
## Warning: Removed 560 rows containing missing values (geom_point).
## Warning: Removed 140 rows containing missing values (geom_point).
## Warning: Removed 756 rows containing missing values (geom_point).
## Warning: Removed 189 rows containing missing values (geom_point).
## Warning: Removed 912 rows containing missing values (geom_point).
## Warning: Removed 228 rows containing missing values (geom_point).
## Warning: Removed 1128 rows containing missing values (geom_point).
## Warning: Removed 282 rows containing missing values (geom_point).
## Warning: Removed 1312 rows containing missing values (geom_point).
## Warning: Removed 328 rows containing missing values (geom_point).
## Warning: Removed 1492 rows containing missing values (geom_point).
## Warning: Removed 373 rows containing missing values (geom_point).
## Warning: Removed 1628 rows containing missing values (geom_point).
## Warning: Removed 407 rows containing missing values (geom_point).
## Warning: Removed 1848 rows containing missing values (geom_point).
## Warning: Removed 462 rows containing missing values (geom_point).
## Warning: Removed 2016 rows containing missing values (geom_point).
## Warning: Removed 504 rows containing missing values (geom_point).
## Warning: Removed 2204 rows containing missing values (geom_point).
## Warning: Removed 551 rows containing missing values (geom_point).
## Warning: Removed 2340 rows containing missing values (geom_point).
## Warning: Removed 585 rows containing missing values (geom_point).
## Warning: Removed 2532 rows containing missing values (geom_point).
## Warning: Removed 633 rows containing missing values (geom_point).
## Warning: Removed 2704 rows containing missing values (geom_point).
## Warning: Removed 676 rows containing missing values (geom_point).
Var = "MSLP"
Varname = "Mean Sea Level Pressure"
Hgt = 0
Fx = 24
for (Fx in Fx_Hour[2:length(Fx_Hour)]) {
fx_subset = CI_Ensemble_Stats %>% filter((Variable == Var) &
(Height == Hgt) &
(Fx_Hour == Fx) ) %>%
select(-Quarter)
seasonal_subset = CI_Ensemble_Stats %>% filter((Variable == Var) &
(Height == Hgt) &
(Fx_Hour == Fx) )
myplot = ggplot(data = seasonal_subset) +
aes(x = Ens_StDev,
y = RMSE_Ens000,
color = Quarter) +
facet_wrap(facets = ~ Quarter) +
theme_bw() +
theme(strip.background = element_rect(fill=NA),
) +
labs(title = str_c(Fx,
"-hr Forecast, CI Seasonal Triangles for ",
Varname,
sep = ""),
subtitle = str_c(region_name,
sep = "")) +
xlab("Ensemble Standard Deviation") +
ylab("Root Mean Squared Error") +
geom_point(data = fx_subset,
color = "grey",
alpha = 0.5) +
scale_colour_manual(values=c("DJF" = "cyan",
"MAM" = "green",
"JJA" = "magenta",
"SON" = "orange"),
guide =FALSE) +
geom_point(data = seasonal_subset,
alpha = 0.7)
print(myplot)
}
## Warning: Removed 340 rows containing missing values (geom_point).
## Warning: Removed 85 rows containing missing values (geom_point).
## Warning: Removed 560 rows containing missing values (geom_point).
## Warning: Removed 140 rows containing missing values (geom_point).
## Warning: Removed 756 rows containing missing values (geom_point).
## Warning: Removed 189 rows containing missing values (geom_point).
## Warning: Removed 912 rows containing missing values (geom_point).
## Warning: Removed 228 rows containing missing values (geom_point).
## Warning: Removed 1128 rows containing missing values (geom_point).
## Warning: Removed 282 rows containing missing values (geom_point).
## Warning: Removed 1312 rows containing missing values (geom_point).
## Warning: Removed 328 rows containing missing values (geom_point).
## Warning: Removed 1492 rows containing missing values (geom_point).
## Warning: Removed 373 rows containing missing values (geom_point).
## Warning: Removed 1628 rows containing missing values (geom_point).
## Warning: Removed 407 rows containing missing values (geom_point).
## Warning: Removed 1852 rows containing missing values (geom_point).
## Warning: Removed 463 rows containing missing values (geom_point).
## Warning: Removed 2016 rows containing missing values (geom_point).
## Warning: Removed 504 rows containing missing values (geom_point).
## Warning: Removed 2204 rows containing missing values (geom_point).
## Warning: Removed 551 rows containing missing values (geom_point).
## Warning: Removed 2340 rows containing missing values (geom_point).
## Warning: Removed 585 rows containing missing values (geom_point).
## Warning: Removed 2536 rows containing missing values (geom_point).
## Warning: Removed 634 rows containing missing values (geom_point).
## Warning: Removed 2704 rows containing missing values (geom_point).
## Warning: Removed 676 rows containing missing values (geom_point).
Var = "M10"
Varname = "10-m Wind Speed"
Hgt = 10
Fx = 24
for (Fx in Fx_Hour[2:length(Fx_Hour)]) {
fx_subset = CI_Ensemble_Stats %>% filter((Variable == Var) &
(Height == Hgt) &
(Fx_Hour == Fx) ) %>%
select(-Quarter)
seasonal_subset = CI_Ensemble_Stats %>% filter((Variable == Var) &
(Height == Hgt) &
(Fx_Hour == Fx) )
myplot = ggplot(data = seasonal_subset) +
aes(x = Ens_StDev,
y = RMSE_Ens000,
color = Quarter) +
facet_wrap(facets = ~ Quarter) +
theme_bw() +
theme(strip.background = element_rect(fill=NA),
) +
labs(title = str_c(Fx,
"-hr Forecast, CI Seasonal Triangles for ",
Varname,
sep = ""),
subtitle = str_c(region_name,
sep = "")) +
xlab("Ensemble Standard Deviation") +
ylab("Root Mean Squared Error") +
geom_point(data = fx_subset,
color = "grey",
alpha = 0.5) +
scale_colour_manual(values=c("DJF" = "cyan",
"MAM" = "green",
"JJA" = "magenta",
"SON" = "orange"),
guide =FALSE) +
geom_point(data = seasonal_subset,
alpha = 0.7)
print(myplot)
}
## Warning: Removed 332 rows containing missing values (geom_point).
## Warning: Removed 83 rows containing missing values (geom_point).
## Warning: Removed 552 rows containing missing values (geom_point).
## Warning: Removed 138 rows containing missing values (geom_point).
## Warning: Removed 748 rows containing missing values (geom_point).
## Warning: Removed 187 rows containing missing values (geom_point).
## Warning: Removed 904 rows containing missing values (geom_point).
## Warning: Removed 226 rows containing missing values (geom_point).
## Warning: Removed 1128 rows containing missing values (geom_point).
## Warning: Removed 282 rows containing missing values (geom_point).
## Warning: Removed 1312 rows containing missing values (geom_point).
## Warning: Removed 328 rows containing missing values (geom_point).
## Warning: Removed 1484 rows containing missing values (geom_point).
## Warning: Removed 371 rows containing missing values (geom_point).
## Warning: Removed 1620 rows containing missing values (geom_point).
## Warning: Removed 405 rows containing missing values (geom_point).
## Warning: Removed 1844 rows containing missing values (geom_point).
## Warning: Removed 461 rows containing missing values (geom_point).
## Warning: Removed 2008 rows containing missing values (geom_point).
## Warning: Removed 502 rows containing missing values (geom_point).
## Warning: Removed 2196 rows containing missing values (geom_point).
## Warning: Removed 549 rows containing missing values (geom_point).
## Warning: Removed 2332 rows containing missing values (geom_point).
## Warning: Removed 583 rows containing missing values (geom_point).
## Warning: Removed 2528 rows containing missing values (geom_point).
## Warning: Removed 632 rows containing missing values (geom_point).
## Warning: Removed 2696 rows containing missing values (geom_point).
## Warning: Removed 674 rows containing missing values (geom_point).
Var = "GUST"
Varname = "10-m Wind Gusts"
Hgt = 10
Fx = 24
for (Fx in Fx_Hour[2:length(Fx_Hour)]) {
fx_subset = CI_Ensemble_Stats %>% filter((Variable == Var) &
(Height == Hgt) &
(Fx_Hour == Fx) ) %>%
select(-Quarter)
seasonal_subset = CI_Ensemble_Stats %>% filter((Variable == Var) &
(Height == Hgt) &
(Fx_Hour == Fx) )
myplot = ggplot(data = seasonal_subset) +
aes(x = Ens_StDev,
y = RMSE_Ens000,
color = Quarter) +
facet_wrap(facets = ~ Quarter) +
theme_bw() +
theme(strip.background = element_rect(fill=NA)) +
labs(title = str_c(Fx,
"-hr Forecast, CI Seasonal Triangles for ",
Varname,
sep = ""),
subtitle = str_c(region_name,
sep = "")) +
xlab("Ensemble Standard Deviation") +
ylab("Root Mean Squared Error") +
geom_point(data = fx_subset,
color = "grey",
alpha = 0.5) +
scale_colour_manual(values=c("DJF" = "cyan",
"MAM" = "green",
"JJA" = "magenta",
"SON" = "orange"),
guide =FALSE) +
geom_point(data = seasonal_subset,
alpha = 0.7)
print(myplot)
}
## Warning: Removed 304 rows containing missing values (geom_point).
## Warning: Removed 76 rows containing missing values (geom_point).
## Warning: Removed 488 rows containing missing values (geom_point).
## Warning: Removed 122 rows containing missing values (geom_point).